Method for rendering on basis of hemispherical orthogonal function

ABSTRACT

The invention discloses a method for rendering on the basis of hemispherical orthogonal function, the method comprising the following steps: selecting rendering fragments and establishing a local coordinate system; acquiring a bidirectional reflectance distribution function of a material; if global illumination is an orthogonal function, determining a rotation matrix of an orthogonal function coefficient according to the rotation angles of the global coordinate system and the local coordinate system, and calculating a local orthogonal function illumination coefficient; converting the local orthogonal function illumination coefficient into a hemispherical orthogonal function illumination coefficient; sampling to obtain the spatial distribution of a bidirectional reflection distribution function of a rendered material; obtaining a hemispherical orthogonal function of the bidirectional reflection distribution function of the rendered material; and using the dot product of a hemispherical orthogonal function coefficient of illumination and a hemispherical orthogonal function coefficient of the bidirectional reflection distribution function of the rendered material and accumulating to obtain the light intensity in the reflection direction. A hemispherical harmonic function(HSH) is used to fit measurement data or theoretically derived data of a BRDF, which may avoid the difficulty of fitting that accurs for a hemispherical harmonics function due to data being missing in the lower hemisphere.

TECHNICAL FIELD

The invention belongs to the technical field of three-dimensional computer image processing, and specifically relates to a rendering method based on a hemispherical orthogonal function.

TECHNICAL BACKGROUND

In computer graphics rendering, the object is discretized into a lot of fragments. According to the material reflection characteristics of the fragments and the incident light emissivity, it is calculated to obtain the reflected emissivity in a specific direction is the core of rendering. This process can be calculated by formula (1):

L _(o)(θ_(o),φ_(o))=∫_(φ) _(i) ₌₀ ^(2π)∫_(θ) _(i) ₌₀ ^(π/2)ƒ_(r)(θ_(i),φ_(i) ^(·),θ_(o),φ_(o))·L _(i)(θ_(i),φ_(i))·cos θ_(i)·sinθ_(i) dθ _(i) dφ _(i)

where:

L_(o) Is along the center P of the element on the object Radiance of reflected light in direction (θ_(o),φ_(o));

L_(i) Is along P Radiance of incident light in direction (θ_(i),φ_(i));

ƒ_(r) Is the bidirectional reflection distribution function BRDF (Bidirectional Reflectance Distribution Function) of the rendered material, which is defined as the ratio of the increment of reflected emissivity dL_(r) (radiance) in a given direction (θ_(r),ϕ_(r)) to the increment of irradiance dE_(i) (irradiance) in the incident direction (θ_(i),ϕ_(i)), see formula (2), which describes the model of how the incident light from different directions is distributed in each exit direction after being reflected by a specific surface;

${BRD{F\left( {\theta_{i},{\phi_{i};\theta_{r}},\phi_{r}} \right)}} = \frac{d{L_{r}\left( {\theta_{i},{\phi_{i};\theta_{r}},\phi_{r}} \right)}}{{dE}_{i}\left( {\theta_{i},\phi_{i}} \right)}$

θ_(i) Is the zenith angle of incident light (in the local coordinate system of the rendered material), The range of θ_(i) is 0 to π/2.

φ_(i) Is the azimuth angle of the incident light (in the local coordinate system of the rendered material), the range of φ_(i) is 0 to 2π.

The integration is performed in the upper half of the rendered fragment, which is a hemispherical space.

The complexity of the reflection characteristics is mainly reflected in the following three aspects:

First of all, BRDF is a function of 4 variables, with high dimensionality, and low efficiency in expression, storage and calculation;

Secondly, it is difficult to calculate the incident light radiance, especially when global illumination is considered. It is necessary to calculate both direct illumination and environmental occlusion and indirect illumination between objects, which requires a large amount of calculation;

Finally, the two-dimensional integration on the hemisphere also requires a large amount of calculation.

Obtaining the reflected emissivity in a specific direction is the core of computer graphics rendering. The currently obtained methods have their own shortcomings in the required time, calculation accuracy, and scope of application.

SUMMARY OF THE INVENTION

The purpose of the implementation of the present invention is to provide a real-time rendering method based on a new hemispherical orthogonal function in view of the shortcomings of the prior art, thereby improving the technical effect of the frame rate and fidelity of rendering.

In order to achieve the above objectives, the present invention provides the following technical solutions:

Method for rendering on the basis of hemispherical orthogonal function, including the following steps:

S1: Select the rendering tablet, establish a local coordinate system;

S2: Get the bidirectional reflection distribution function of the material;

S3: Determine if the global light is an orthogonal function; if it is, step S4 and S5; if not, proceeds to step S4′ and S5′;

S4: Determine the rotation matrix of the orthogonal function coefficient according to the corner of the global coordinate system and the local coordinate system, and calculate the local orthogonal function of the local orthogonal function;

S5: According to the equivalent formula, the local orthogonal function is converted to a hemispheric orthogonal function; the hemisphere orthogonal function includes three types of hemispherical harmonic functions, ZerniKe hemisphere functions, and vertical hemisphere functions;

S4′: Sampling of the upper half of the ball, obtaining the spatial distribution of light;

S5′: Fit light distribution with hemisphere orthogonal functions, obtains the slime orthogonal function of the light; the hemispicular orthogonal function includes a hemispherical harmonic function, the zernike hemisphere function, and the vertical hemisphere function;

S6: Determines the direction of the reflected light, the spatial distribution of the two-way reflection distribution function of the rendered material is obtained;

S7: Fit the two-way reflection distribution function of the rendered material with a hemisphere orthogonal function, and obtain a hemisphere orthogonal function factor of the two-way reflection distribution function of the rendered material;

S8: The hemisphere orthogonal function coefficient of light is multiplied by the half-way reflection distribution function of the rendered material and accumulates the light intensity of the reflection direction.

Furthermore, the method is characterized in that With a hemispheric orthogonal function, the bidirectional distribution function of the incident light and the material is converted from the airspace to the frequency domain, and the radiation in the reflection direction is obtained by multiplication in the frequency domain.

Furthermore, the method is characterized in that The hemispherical harmonic function HSH:

${{HS{H_{n}^{m}\left( {s,\varphi} \right)}} = {\left( {- 1} \right)^{❘m❘}{\sqrt{\frac{\left( {{2n} + 1} \right)}{2}\frac{\left( {n - {❘m❘}} \right)!}{\left( {n + {❘m❘}} \right)!}} \cdot P_{n}^{❘m❘}}{\left\{ {\cos(\theta)} \right\} \cdot {\Phi_{m}(\varphi)}}}},$ $\begin{pmatrix} {{n = {0,2,4}},{\ldots;{m = {{0,} \pm 2}}},\ldots,{\pm n}} \\ {{n = {1,3,5}},{\ldots;{m = {{{\pm 1},} \pm 3}}},\ldots,{\pm n}} \end{pmatrix}$

Where, θ is the top angle; φ Is a square angle;

${P_{n}^{❘m❘}(x)} = {\frac{\left( {- 1} \right)^{❘m❘}}{2^{n}{n!}}\left( {1 - x^{2}} \right)^{\frac{❘m❘}{2}}\left( \frac{d}{dx} \right)^{n + {❘m❘}}{\left( {x^{2} - 1} \right)^{n}.}}$

Furthermore, the method is characterized in that The Zernik Hemisphere Function ZSF:

${{ZS{F_{n}^{m}\left( {t,\varphi} \right)}} = {\sqrt{2\left( {n + 1} \right)} \cdot {R_{n}^{m}(t)} \cdot {\Phi_{m}(\varphi)}}},$ $\begin{pmatrix} {{n = {0,2,4}},{\ldots;{m = {{0,} \pm 2}}},\ldots,{\pm n}} \\ {{n = {1,3,5}},{\ldots;{m = {{{\pm 1},} \pm 3}}},\ldots,{\pm n}} \end{pmatrix}$ ${Where},{{R_{n}^{m}(t)} = {\sum\limits_{k = 0}^{{({n - {❘m❘}})}/2}{\frac{\left( {- 1} \right)^{k}{\left( {n - k} \right)!}}{{k!}{\left( {\frac{n + {❘m❘}}{2} - k} \right)!}{\left( {\frac{n - {❘m❘}}{2} - k} \right)!}}{t^{n - {2k}}.}}}}$

Furthermore, the method is characterized in that The longitudinal hemisphere function LSF:

${{LS{F_{n}^{m}\left( {w,\varphi} \right)}} = {\sqrt{\frac{2^{{❘m❘} + 2.5}}{\gamma_{v}^{({a,b})}}}w^{❘m❘}{{J\left( {\frac{n - {❘m❘}}{2},0,\frac{{2{❘m❘}} - 1}{2},{{2w^{2}} - 1}} \right)} \cdot {\Phi_{m}(\varphi)}}}},$ $\begin{pmatrix} {{n = {0,2}},{\ldots;{m = {{0,} \pm 2}}},\ldots,{\pm n}} \\ {{n = {1,3}},{\ldots;{m = {{{\pm 1},} \pm 3}}},\ldots,{\pm n}} \end{pmatrix}$ ${Where},{{J\left( {v,a,b,x} \right)} = {\frac{1}{v!}{\sum\limits_{k = 0}^{v}{\frac{\left( {- v} \right)_{k}\left( {a + b + v + 1} \right)_{k}\left( {a + k + 1} \right)_{v - k}}{k!}{\left( \frac{1 - x}{2} \right)^{k}.}}}}}$

Compared with the prior art, the present invention has the following advantages:

1. The expression of BRDF and incident light is more effective, and the fitting process is numerically stable. Using the hemispherical harmonic function HSH to fit BRDF measurement data or theoretically derived data can avoid the difficulty of fitting the spherical harmonic function due to the missing data in the lower hemisphere (the fitting result has a large amplitude and is sensitive to measurement noise).

2. If global illumination is expressed by spherical harmonic function, this method provides a fast equivalent conversion method of global illumination spherical harmonic function coefficients to hemispherical illumination function coefficients in the local coordinate system, which improves the calculation speed.

3. The Fourier transform on the hemispherical surface is used to calculate the reflected illuminance in the frequency domain. Compared with the method of convolution in the spatial domain, the multiplication in the frequency domain has the advantages of less calculation, fast speed, and the cutoff frequency can be selected according to the fineness of rendering.

DESCRIPTION OF THE FIGURES

FIG. 1 illustrates a flowchart of image rendering algorithm based on hemispherical functions.

FIG. 2 illustrates the definition of bidirectional reflectance distribution function (BRDF).

FIG. 3 is a flowchart of calculation of surface element reflection light based on hemispherical orthogonal function.

DETAILED DESCRIPTION

The rendering method based on the hemispherical orthogonal function of the present invention mainly includes:

(1) Select the patch for rendering, obtain its position, orientation and material reflectance characteristics. (2) Using the hemispherical functions to fit the incident light and obtain its representation by coefficients (3) Using the hemispherical functions to fit the bidirectional distribution function along specific perspective direction and obtain its representation by coefficients. (4) Based on the orthogonal property of hemispherical functions, the two-dimensional integral in rendering equation can be calculated by dot product of illumination coefficients and BRDF coefficients.

Here the hemispherical functions are hemispherical harmonics (HSH), Zernkie hemispherical functions (ZHF) and Longitudinal hemispherical functions (LHF).

Optimizely, in step (2) if the illumination is not expressed by spherical harmonics, operations in step (2-2) will be carried. Otherwise, operations in step (2-2) will be carried.

(2-1) Ray tracing is used to calculate the intensity of illumination light from different directions to the rendering patch over its upper hemisphere. Then the distribution is fitted by hemispherical functions and coefficients of illumination are obtained.

(2-2) The illumination has been expressed as spherical functions. HSH is preferred to use. An approach is present to efficiently convert spherical function coefficients in global coordinates to HSH coefficients in local coordinate system. The specific implementation of the present invention will be further described below in conjunction with the drawings and embodiments.

The detailed process to perform the presented method will be described using diagrams. FIG. 1 illustrates a flowchart of image rendering algorithm based on hemispherical functions.

Step 1: Select the patch for rendering, establish local coordinates systems for rendering.

Step 1.1: Obtain patch center position P in space {x_(p),y_(p),z_(p)}, its normal direction {right arrow over (n)} and material BRDF. The definition of BRDF is illustrated in FIG. 2 where dA is patch being rendered, {right arrow over (N)} is the normal direction of the patch, {right arrow over (x)} is the positive x direction for local coordinates. {right arrow over (y)} is the positive y direction for local coordinates. {right arrow over (z)} is the positive z direction for local coordinates. θ_(i) is the polar angle of incident light. φ_(i) is the azimuth angle of incident light. dω_(i) is the differential solid angle along incident direction. θ_(r) is the polar angle of reflectance light. φ_(r) is the azimuth angle of reflectance light. dω_(r) is the differential solid angle along reflectance direction.

Step 1.2: Establish the local rendering coordinates system. Position P is the origin. Normal direction {right arrow over (n)} is the Z positive direction. Two perpendicular direction X and Y are selected to match with BRDF directions to simplify the computation.

Step 2: Using the hemispherical functions to fit the incident light and obtain its representation by coefficients. Hemispherical functions are hemispherical harmonics (HSH), Zernkie hemispherical functions (ZHF) and Longitudinal hemispherical functions (LHF), in which:

(1) The expression of the hemispherical harmonic function (HSH) is as follows:

${{HS{H_{n}^{m}\left( {s,\varphi} \right)}} = {\left( {- 1} \right)^{❘m❘}{\sqrt{2\left( {{2n} + 1} \right)\frac{\left( {n - {❘m❘}} \right)!}{\left( {n + {❘m❘}} \right)!}} \cdot P_{n}^{❘m❘}}{\left\{ {\cos(\theta)} \right\} \cdot {\Phi_{m}(\varphi)}}}},$ $\begin{pmatrix} {{n = {0,2,4}},{\ldots;{m = {{0,} \pm 2}}},\ldots,{\pm n}} \\ {{n = {1,3,5}},{\ldots;{m = {{{\pm 1},} \pm 3}}},\ldots,{\pm n}} \end{pmatrix}$

among them P_(n) ^(|m|) Is associated Legendre Polynomials (Associated Legendre

Polynomials), which is expressed as follows:

${P_{n}^{❘m❘}(x)} = {\frac{1}{2^{n}{n!}}\left( {1 - x^{2}} \right)^{\frac{❘m❘}{2}}\left( \frac{d}{dx} \right)^{n + {❘m❘}}\left( {x^{2} - 1} \right)^{n}}$

(2) The expression of Zernkie hemispherical functions (ZHF) is as follows:

${{{ZSF}_{n}^{m}\left( {t,\varphi} \right)} = {\sqrt{2\left( {n + 1} \right)} \cdot {R_{n}^{m}(t)} \cdot {\Phi_{m}(\varphi)}}},\begin{pmatrix} {{n = 0},2,4,{\cdots;{m = 0}},{\pm 2},\cdots,{\pm n}} \\ {{n = 1},3,5,{{\cdots;} = {\pm 1}},{\pm 3},\cdots,{\pm n}} \end{pmatrix}$

among them R_(n) ^(m)(t) Is a Zernike polynomial, expressed as follows:

${{{R_{n}^{m}(t)} = {\overset{{({n - {❘m❘}})}/2}{\sum\limits_{k = 0}}{\frac{\left( {- 1} \right)^{k}{\left( {n - k} \right)!}}{{k!}{\left( {\frac{n + {❘m❘}}{2} - k} \right)!}{\left( {\frac{n - {❘m❘}}{2} - k} \right)!}}t^{n - {2k}}}}};}{{w = {1 - {\cos(\theta)}}};}$

(3) The expression of the longitudinal spherical function (LSF) is as follows:

${{{LSF}_{n}^{m}\left( {w,\varphi} \right)} = {\sqrt{\frac{2^{{❘m❘} + 2.5}}{\gamma_{v}^{({a,b})}}}w^{❘m❘}{{J\left( {\frac{n - {❘m❘}}{2},0,\frac{{2{❘m❘}} - 1}{2},{{2w^{2}} - 1}} \right)} \cdot {\Phi_{m}(\varphi)}}}},\begin{pmatrix} {{n = 0},2,{\cdots;{m = 0}},{\pm 2},\cdots,{\pm n}} \\ {{n = 1},3,{\cdots;{m = {\pm 1}}},{\pm 3},\cdots,{\pm n}} \end{pmatrix}$

Where J is the Jacobi polynomial, as shown below:

${J\left( {v,a,b,x} \right)} = {\frac{1}{v!}{\overset{v}{\sum\limits_{k = 0}}{\frac{\left( {- v} \right)_{k}\left( {a + b + v + 1} \right)_{k}\left( {a + k + 1} \right)_{v - k}}{k!}\left( \frac{1 - x}{2} \right)^{k}}}}$

This embodiment only uses the hemispherical harmonic function HSH as an example to describe the method of the present invention in detail. The application methods of the Zernike hemispherical function ZSF and the longitudinal hemispherical function LSF are the same as the hemispherical harmonic function HSH, and will not be repeated here.

Step 2.1: Calculate the illumination HSH coefficients by fitting

Step 2.1.1: Discrete the upper hemispherical space into different directions. The discrete number is s. Use ray tracing or voxels tracing approach, the illumination intensities along these directions are obtained. Use a vector to record the result as:

{right arrow over (L)}_(i)=[l_(i) ¹ l_(i) ² ⋅⋅⋅ l_(i) ^(s)]^(T)

Step 2.1.2: Use HSH functions degree no more than n to fit the illumination light and obtain the fitting coefficients. The total number of functions is N (N=(n+1)*(n+2)/2). Use a vector to record the results.

{right arrow over (A)}=[a₁ a₂ ⋅⋅⋅ a_(N)]hu T

And the incident light can be approximated as following.

${L_{i}\left( {\theta_{i},\varphi_{i}} \right)} = {\overset{N}{\sum\limits_{l = 1}}{a_{l} \cdot {{HSH}_{l}\left( {\theta_{i},\varphi_{i}} \right)}}}$

Step 2.2: Calculate the illumination HSH coefficients by rotating SH. In some circumstance, the illumination is represented by spherical harmonics. An approach is put forward to make an efficient transfer.

Step 2.2.1: Rotate the global illumination to local rendering coordinates system. The global illumination is recorded as a set of spherical harmonic coefficients as following.

{right arrow over (E)}=[ε₀ ¹ ε₀ ² ⋅⋅⋅ ε₀ ^(N)]^(T)

The Euler angle between global system and local rendering system can be calculated and denoted as α, β, γ. Let name the illumination coefficients in local coordinates as {right arrow over (E)}_(r).

{right arrow over (E)}_(r)=[ε_(r) ¹ ε_(r) ² ⋅⋅⋅ ε_(r) ^(N)]^(T)

Then the transfer calculation can be carried by equation.

{right arrow over (E)} _(r) =R _(SH)(α,β,γ)·{right arrow over (E)} ₀

where R_(SH)(α,β,γ) is the rotational matrix for spherical harmonics. It can be derived by recursive method, or by ZXZXZ approach.

Step 2.2.2: Transfer the local illumination SH coefficients into equivalent HSH coefficients. The formula of equivalent transfer is listed.

$\left\{ \begin{matrix} {{{SH}_{n}^{m}\left( {\theta,\varphi} \right)} = {\overset{\infty}{\sum\limits_{k = 0}}{\gamma_{n,m}^{{m + {2k}},m} \cdot {{HSH}_{m + {2k}}^{m}\left( {\theta,\varphi} \right)}}}} \\ {\gamma_{n,m}^{{m + {2k}},m} = {\left\langle {{SH}_{n}^{m},{HSH}_{m + {2k}}^{m}} \right\rangle/\left( {2\pi} \right)}} \end{matrix} \right.$

Where

U,V

is the integral of function of U and V over hemispherical surface.

U,V

≡∫₀ ^(π/2)∫₀ ^(2π) U(θ,φ)·V(θ,φ)·sinθ·dθ·dφ

If (n−m) is an even number, SH_(n) ^(m)(θ,φ)=HSH_(n) ^(m)(θ,φ), the coefficient of SH_(n) ^(m)(θ,φ) equals to the coefficient of HSH_(n) ^(m)(θ,φ). If (n−m) is odd, the calculation of transferring factor γ_(n,m) ^(m+2k,m) is carried to perform the equivalent transfer calculation.

For example,

${{SH}_{1}^{0}\left( {\theta,\varphi} \right)} = {{\frac{\sqrt{3}}{2}{{HSH}_{0}^{0}\left( {\theta,\varphi} \right)}} + {\frac{\sqrt{15}}{8}{{HSH}_{2}^{0}\left( {\theta,\varphi} \right)}} - {\frac{\sqrt{3}}{16}{{HSH}_{4}^{0}\left( {\theta,\varphi} \right)}} + {\frac{\sqrt{39}}{128}{{HSH}_{6}^{0}\left( {\theta,\varphi} \right)}} - \cdots}$ ${{SH}_{2}^{1}\left( {\theta,\varphi} \right)} = {{\frac{3\sqrt{5}}{8}{{HSH}_{1}^{1}\left( {\theta,\varphi} \right)}} + {\frac{\sqrt{70}}{16}{{HSH}_{3}^{1}\left( {\theta,\varphi} \right)}} - {\frac{5\sqrt{11}}{128}{{HSH}_{5}^{1}\left( {\theta,\varphi} \right)}} + {\frac{\sqrt{63}}{128}{{HSH}_{7}^{1}\left( {\theta,\varphi} \right)}} - \cdots}$

Set an truncation to get an acceptable transfer precision.

Step 2.2.3: Combining like items, equivalent HSH coefficients are Calculated.

$a_{n}^{m} = {\overset{\infty}{\sum\limits_{k = 0}}{\varepsilon_{r}^{{m + k},m} \cdot \gamma_{{m + k},m}^{n,m}}}$

among them ε_(r) ^(m+k,m) Is the coefficient of spherical harmonic illumination function in the local coordinate system, take an example as coefficient a₀ ⁰ for HSH₀ ⁰(θ,φ)

${a_{1} = {a_{0}^{0} = {{\varepsilon_{r}^{0,0} \cdot \gamma_{0,0}^{0,0}} + {\varepsilon_{r}^{1,0} \cdot \gamma_{1,0}^{0,0}} + {\varepsilon_{r}^{2,0} \cdot \gamma_{2,0}^{0,0}} + {\varepsilon_{r}^{3,0} \cdot \gamma_{3,0}^{0,0}} + \cdots}}}{= {{1. \cdot \varepsilon_{r}^{0,0}} + {\frac{\sqrt{3}}{2} \cdot \varepsilon_{r}^{1,0}} + {0. \cdot \varepsilon_{r}^{2,0}} - {\frac{\sqrt{7}}{8} \cdot \varepsilon_{r}^{3,0}} + \cdots}}$

Optimization1: Since HSH and SH have the same part in expression, half of the conversion coefficient is a one-to-one conversion, which is beneficial to calculation. ZSF and LSF do not have this advantage, and the calculation speed will be slower, so HSH is preferred.

Optimization2: The transfer matrix from SH to HSH is fixed, and the conversion coefficient is calculated in advance and does not need to be performed in real time. Therefore, this conversion requires less calculation and fast execution speed. So far, regardless of the calculation method 1 or 2, the hemispheric function fitting coefficient of the incident light irradiance {right arrow over (A)}=[a₁ a₂ ⋅⋅⋅a_(N)]^(T) can be obtained.

Step 3: Using the hemispherical functions to fit the bidirectional distribution function along specific perspective direction and obtain its representation by coefficients. The bidirectional reflection distribution function (BRDF) is defined as the target in a certain direction (θ_(r),ϕ_(r)) On the reflection brightness dL_(r) and incident direction (θ_(i),ϕ_(i)) Illuminance dE_(i) ratio. The BRDF of the material can be calculated by surface physical modeling, and the actual measurement can also be carried out by instruments.

Step 3.1: Determine the expression in the local coordinate system of the direction to calculate the reflected emissivity (θ_(r),φ_(r))

Step 3.2: Fixing (θ_(r),φ_(r)), then through physical models or experimental data, the influence of irradiance in each incident direction (θ_(i),φ_(i)) on the emissivity in the reflection direction (θ_(r), φ_(r)) is obtained, recorded as BRDF(θ_(i),φ_(i))·cos(θ_(i))|(θ_(r),φ_(r)). Obviously BRDF(θ_(i),φ_(i))·(θ_(i)) Is a function of the hemispherical distribution. If it has t sampling directions in the hemispherical space, the expression in the spatial domain is {right arrow over (F)}_(r)=[ƒ_(r) ¹ ƒ_(r) ² ⋅⋅⋅ ƒ_(r) ^(t)]^(T).

Step 3.3: Fit {right arrow over (F)}_(r) with HSH functions degree no more than n (N terms), and obtain fitting coefficient {right arrow over (B)}=[b₁ b₂ ⋅⋅⋅ b_(N)]^(T), BRDF can use the first N HSH to express:

${{{BRDF}\left( {\theta_{i},\varphi_{i}} \right)}{\cos\left( \theta_{i} \right)}} = {\overset{N}{\sum\limits_{l = 1}}{b_{l} \cdot {{HSH}_{l}\left( {\theta_{i},\varphi_{i}} \right)}}}$

Step 4: Substitute the reflection integral formula to calculate Reflected lighting of the direction (θ_(r),φ_(r)), the calculation process of surface element reflected light based on hemispherical orthogonal function is shown in FIG. 3 , The specific method is as follows:

L_(o)(θ_(o), φ_(o))  = ∫_(φ_(i) = 0)^(2π)∫_(θ_(i) = 0)^(π/2)f_(r)(θ_(i), φ_(i); θ_(o), φ_(o)) ⋅ L_(i)(θ_(i), φ_(i)) ⋅ cos θ_(i) ⋅ sin θ_(i)dθ_(i)dφ_(i) $= {\int_{\varphi_{i} = 0}^{2\pi}{\int_{\theta_{i} = 0}^{\pi/2}\left( {\overset{N}{\sum\limits_{k = 1}}{{a_{k} \cdot {{HSH}_{k}\left( {\theta_{i},\varphi_{i}} \right)} \cdot \left( {\overset{N}{\sum\limits_{l = 1}}{b_{l} \cdot {{HSH}_{l}\left( {\theta_{i},\varphi_{i}} \right)}}} \right) \cdot \sin}\theta_{i}d\theta_{i}d\varphi_{i}}} \right.}}$ $\left. {{\left. {= {\overset{N}{\sum\limits_{k = 1}}{\overset{N}{\sum\limits_{l = 1}}{\left( {a_{k}b_{l}} \right) \cdot \left\{ {\int_{\varphi_{i} = 0}^{2\pi}{\int_{\theta_{i} = 0}^{\pi/2}{{{HSH}_{k}\left( {\theta_{i},\varphi_{i}} \right)} \cdot {{HSH}_{l}\left( {\theta_{i},\varphi_{i}} \right)}}}} \right.}}}} \right) \cdot \sin}\theta_{i}d\theta_{i}d\varphi_{i}} \right\}$ $= {2\pi{\overset{N}{\sum\limits_{l = 1}}{a_{l}b_{l}}}}$

The above descriptions are only preferred embodiments of the present invention, and are not intended to limit the present invention. Any modification, equivalent replacement and improvement made within the spirit and principle of the present invention shall be included in the protection scope of the present invention. 

1. A method for rendering on the basis of hemispherical orthogonal function, including the following steps: S1: Select the rendering tablet, establish a local coordinate system; S2: Get the bidirectional reflection distribution function of the material; S3: Determine if the global light is an orthogonal function; if it is, step S4 and S5; if not, proceeds to step S4′ and S5′; S4: Determine the rotation matrix of the orthogonal function coefficient according to the corner of the global coordinate system and the local coordinate system, and calculate the local orthogonal function of the local orthogonal function; S5: According to the equivalent formula, the local orthogonal function is converted to a hemispheric orthogonal function; the hemisphere orthogonal function includes three types of hemispherical harmonic functions, ZerniKe hemisphere functions, and vertical hemisphere functions; S4′: Sampling of the upper half of the ball, obtaining the spatial distribution of light; S5′: Fit light distribution with hemisphere orthogonal functions, obtains the slime orthogonal function of the light; the hemispicular orthogonal function includes a hemispherical harmonic function, the zernike hemisphere function, and the vertical hemisphere function; S6: Determines the direction of the reflected light, the spatial distribution of the two-way reflection distribution function of the rendered material is obtained; S7: Fit the two-way reflection distribution function of the rendered material with a hemisphere orthogonal function, and obtain a hemisphere orthogonal function factor of the two-way reflection distribution function of the rendered material; S8: The hemisphere orthogonal function coefficient of light is multiplied by the half-way reflection distribution function of the rendered material and accumulates the light intensity of the reflection direction.
 2. The method for rendering on the basis of hemispherical orthogonal function according to claim 1, wherein with a hemispheric orthogonal function, the bidirectional distribution function of the incident light and the material is converted from the airspace to the frequency domain, and the radiation in the reflection direction is obtained by multiplication in the frequency domain.
 3. The method for rendering on the basis of hemispherical orthogonal function according to claim 1, wherein the hemispherical harmonic function HSH: ${{{HSH}_{n}^{m}\left( {s,\varphi} \right)} = {\left( {- 1} \right)^{❘m❘}{\sqrt{\frac{\left( {{2n} + 1} \right)}{2}\frac{\left( {n - {❘m❘}} \right)!}{\left( {n + {❘m❘}} \right)!}} \cdot P_{n}^{❘m❘}}{\left\{ {\cos(\theta)} \right\} \cdot {\Phi_{m}(\varphi)}}}},\begin{pmatrix} {{n = 0},2,4,{\cdots;{m = 0}},{\pm 2},\cdots,{\pm n}} \\ {{n = 1},3,5,{\cdots;{m = {\pm 1}}},{\pm 3},\cdots,{\pm n}} \end{pmatrix}$ Where, θ is the top angle; φ Is a square angle; ${P_{n}^{❘m❘}(x)} = {\frac{\left( {- 1} \right)^{❘m❘}}{2^{n}{n!}}\left( {1 - x^{2}} \right)^{\frac{❘m❘}{2}}\left( \frac{d}{dx} \right)^{n + {❘m❘}}{\left( {x^{2} - 1} \right)^{n}.}}$
 4. The method for rendering on the basis of hemispherical orthogonal function according to claim 1, wherein the Zernik Hemisphere Function ZSF: ${{{ZSF}_{n}^{m}\left( {t,\varphi} \right)} = {\sqrt{2\left( {n + 1} \right)} \cdot {R_{n}^{m}(t)} \cdot {\Phi_{m}(\varphi)}}},\begin{pmatrix} {{n = 0},2,4,{\cdots;{m = 0}},{\pm 2},\cdots,{\pm n}} \\ {{n = 1},3,5,{{\cdots;} = {\pm 1}},{\pm 3},\cdots,{\pm n}} \end{pmatrix}$ ${R_{n}^{m}(t)} = {\overset{{({n - {❘m❘}})}/2}{\sum\limits_{k = 0}}{\frac{\left( {- 1} \right)^{k}{\left( {n - k} \right)!}}{{k!}{\left( {\frac{n + {❘m❘}}{2} - k} \right)!}{\left( {\frac{n - {❘m❘}}{2} - k} \right)!}}t^{n - {2k}}}}$ Where, .
 5. The method for rendering on the basis of hemispherical orthogonal function according to claim 1, wherein the longitudinal hemisphere function LSF: ${{{LSF}_{n}^{m}\left( {w,\varphi} \right)} = {\sqrt{\frac{2^{{❘m❘} + 2.5}}{\gamma_{v}^{({a,b})}}}w^{❘m❘}{{J\left( {\frac{n - {❘m❘}}{2},0,\frac{{2{❘m❘}} - 1}{2},{{2w^{2}} - 1}} \right)} \cdot {\Phi_{m}(\varphi)}}}},\begin{pmatrix} {{n = 0},2,{\cdots;{m = 0}},{\pm 2},\cdots,{\pm n}} \\ {{n = 1},3,{\cdots;{m = {\pm 1}}},{\pm 3},\cdots,{\pm n}} \end{pmatrix}$ ${J\left( {v,a,b,x} \right)} = {\frac{1}{v!}{\overset{v}{\sum\limits_{k = 0}}{\frac{\left( {- v} \right)_{k}\left( {a + b + v + 1} \right)_{k}\left( {a + k + 1} \right)_{v - k}}{k!}\left( \frac{1 - x}{2} \right)^{k}}}}$ Where, 